2.5 Young/old gene expression
Mean expression of young genes.
We check whether newer genes are restricted in expression during the ‘hourglass’ stage.
library(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.3 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.4 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.0
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(myTAI)Load PES
Non-transformed
Fd_PES <-
readr::read_csv(file = "data/Fd_PES.csv") %>%
dplyr::select(-c("gamete", "matSP"))## Rows: 7907 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (9): PS, gamete, E1, E2, E3, E4, E5, E6, matSP
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.Fs_PES <-
readr::read_csv(file = "data/Fs_PES.csv") %>%
dplyr::select(-c("gamete", "matSP"))## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): PS, gamete, 24H, 48H, 1w, 3w, 4w, matSP
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.sqrt-tranformed
Fd_PES.sqrt <-
readr::read_csv(file = "data/Fd_PES.sqrt.csv") %>%
dplyr::select(-c("gamete", "matSP"))## Rows: 7907 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (9): PS, gamete, E1, E2, E3, E4, E5, E6, matSP
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.Fs_PES.sqrt <-
readr::read_csv(file = "data/Fs_PES.sqrt.csv") %>%
dplyr::select(-c("gamete", "matSP"))## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): PS, gamete, 24H, 48H, 1w, 3w, 4w, matSP
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.log2-tranformed
Fd_PES.log2 <-
readr::read_csv(file = "data/Fd_PES.log2.csv") %>%
dplyr::select(-c("gamete", "matSP"))## Rows: 7907 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (9): PS, gamete, E1, E2, E3, E4, E5, E6, matSP
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.Fs_PES.log2 <-
readr::read_csv(file = "data/Fs_PES.log2.csv") %>%
dplyr::select(-c("gamete", "matSP"))## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): PS, gamete, 24H, 48H, 1w, 3w, 4w, matSP
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.rank-tranformed
Fd_PES.rank <-
readr::read_csv(file = "data/Fd_PES.rank.csv") %>%
dplyr::select(-c("gamete", "matSP"))## Rows: 7907 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (9): PS, gamete, E1, E2, E3, E4, E5, E6, matSP
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.Fs_PES.rank <-
readr::read_csv(file = "data/Fs_PES.rank.csv") %>%
dplyr::select(-c("gamete", "matSP"))## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): PS, gamete, 24H, 48H, 1w, 3w, 4w, matSP
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.rlog-tranformed
Fd_PES.rlog <-
readr::read_csv(file = "data/Fd_PES.rlog.csv") %>%
dplyr::select(-c("gamete", "matSP"))## Rows: 7907 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (9): PS, gamete, E1, E2, E3, E4, E5, E6, matSP
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.Fs_PES.rlog <-
readr::read_csv(file = "data/Fs_PES.rlog.csv") %>%
dplyr::select(-c("gamete", "matSP"))## Rows: 8291 Columns: 9
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (8): PS, gamete, 24H, 48H, 1w, 3w, 4w, matSP
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.Average expression of young genes
I now plot the mean expression of each PS.
plot_PS_expression <- function(PES, stage_names, average = "mean"){
# PES_mean_expression <-
# PES %>%
# tidyr::pivot_longer(!c(PS, GeneID), values_to = "Expression", names_to = "Stage") %>%
# dplyr::group_by(Stage) %>%
# dplyr::summarise(stage_mean_exp = mean(Expression))
averaging_function <- base::match.fun(average)
PES_in <-
PES %>%
tidyr::pivot_longer(!c(PS, GeneID), values_to = "Expression", names_to = "Stage") %>%
# dplyr::left_join(PES_mean_expression) %>%
dplyr::mutate(PS_name = case_when(
PS <=6 ~ "Old",
PS == 7 ~ "Brown Algal",
PS == 8 ~ "Species specific"
)) %>%
dplyr::filter(!PS_name == "Old") %>%
dplyr::group_by(PS_name,Stage) %>%
dplyr::summarise(mean_expression = averaging_function(Expression))
PES_in$Stage <- factor(PES_in$Stage, levels = unique(stage_names))
PES_in$PS_name <- factor(PES_in$PS_name, levels = c("Old", "Brown Algal", "Species specific"))
return(PES_in)
}Fd_stage_names = colnames(Fd_PES)[-c(1,2)]
Fd_PES.sqrt %>%
plot_PS_expression(stage_names = Fd_stage_names) %>%
ggplot2::ggplot(aes(x = Stage, y = mean_expression, group = PS_name, colour = as.factor(PS_name))) +
ggplot2::geom_line(linewidth = 3) +
ggplot2::theme(legend.position = "bottom") +
ggplot2::labs(title = "F. distichus", y = "mean sqrt(TPM)", colour = "gene age") +
ggsci::scale_colour_npg()## `summarise()` has grouped output by 'PS_name'. You can override using the
## `.groups` argument.
Fd_PES.log2 %>%
plot_PS_expression(stage_names = Fd_stage_names) %>%
ggplot2::ggplot(aes(x = Stage, y = mean_expression, group = PS_name, colour = as.factor(PS_name))) +
ggplot2::geom_line(linewidth = 3) +
ggplot2::theme(legend.position = "bottom") +
ggplot2::labs(title = "F. distichus", y = "mean log(TPM)", colour = "gene age") +
ggsci::scale_colour_npg()## `summarise()` has grouped output by 'PS_name'. You can override using the
## `.groups` argument.
Fd_PES.sqrt %>%
plot_PS_expression(stage_names = Fd_stage_names, average = "median") %>%
ggplot2::ggplot(aes(x = Stage, y = mean_expression, group = PS_name, colour = as.factor(PS_name))) +
ggplot2::geom_line(linewidth = 3) +
ggplot2::theme(legend.position = "bottom") +
ggplot2::labs(title = "F. distichus", y = "median sqrt(TPM)", colour = "gene age") +
ggsci::scale_colour_npg()## `summarise()` has grouped output by 'PS_name'. You can override using the
## `.groups` argument.
Fd_PES.log2 %>%
plot_PS_expression(stage_names = Fd_stage_names, average = "median") %>%
ggplot2::ggplot(aes(x = Stage, y = mean_expression, group = PS_name, colour = as.factor(PS_name))) +
ggplot2::geom_line(linewidth = 3) +
ggplot2::theme(legend.position = "bottom") +
ggplot2::labs(title = "F. distichus", y = "median log(TPM)", colour = "gene age") +
ggsci::scale_colour_npg()## `summarise()` has grouped output by 'PS_name'. You can override using the
## `.groups` argument.
Fs_stage_names = colnames(Fs_PES)[-c(1,2)]
Fs_PES.sqrt %>%
plot_PS_expression(stage_names = Fs_stage_names) %>%
ggplot2::ggplot(aes(x = Stage, y = mean_expression, group = PS_name, colour = as.factor(PS_name))) +
ggplot2::geom_line(linewidth = 3) +
ggplot2::theme(legend.position = "bottom") +
ggplot2::labs(title = "F. serratus", y = "mean sqrt(TPM)", colour = "gene age") +
ggsci::scale_colour_npg()## `summarise()` has grouped output by 'PS_name'. You can override using the
## `.groups` argument.
Fs_PES.log2 %>%
plot_PS_expression(stage_names = Fs_stage_names) %>%
ggplot2::ggplot(aes(x = Stage, y = mean_expression, group = PS_name, colour = as.factor(PS_name))) +
ggplot2::geom_line(linewidth = 3) +
ggplot2::theme(legend.position = "bottom") +
ggplot2::labs(title = "F. serratus", y = "mean log(TPM)", colour = "gene age") +
ggsci::scale_colour_npg()## `summarise()` has grouped output by 'PS_name'. You can override using the
## `.groups` argument.
Fs_PES.sqrt %>%
plot_PS_expression(stage_names = Fs_stage_names, average = "median") %>%
ggplot2::ggplot(aes(x = Stage, y = mean_expression, group = PS_name, colour = as.factor(PS_name))) +
ggplot2::geom_line(linewidth = 3) +
ggplot2::theme(legend.position = "bottom") +
ggplot2::labs(title = "F. serratus", y = "median sqrt(TPM)", colour = "gene age") +
ggsci::scale_colour_npg()## `summarise()` has grouped output by 'PS_name'. You can override using the
## `.groups` argument.
Fs_PES.log2 %>%
plot_PS_expression(stage_names = Fs_stage_names, average = "median") %>%
ggplot2::ggplot(aes(x = Stage, y = mean_expression, group = PS_name, colour = as.factor(PS_name))) +
ggplot2::geom_line(linewidth = 3) +
ggplot2::theme(legend.position = "bottom") +
ggplot2::labs(title = "F. serratus", y = "median log(TPM)", colour = "gene age") +
ggsci::scale_colour_npg()## `summarise()` has grouped output by 'PS_name'. You can override using the
## `.groups` argument.
Average expression of all genes
To have an understanding of the variance in the median calculation, I will obtain the standard deviation via bootstrapping. Bootstrapping is done because the distribution of the standard deviation is non-normal.
Turns out it would take longer than I thought.
# Average.preplot <- function(ExpressionSet, permutations = 500, is_PES = F, boot_function = "median"){
# # Check if the input is correct
# print("colnames should be either (1) GeneID and conditions, (2) conditions with GeneIDs as rownames, (3) PS, GeneID and conditions")
#
# if(!boot_function %in% c("mean", "median")){
# message("boot_function must be either 'mean' or 'median' (for now)")
# }
#
# if(is_PES){
# ExpressionSetInput <- dplyr::select(ExpressionSet, -1)
# }
# else{
# if("PS" %in% colnames(ExpressionSet)){
# stop("Check the is_PES option or if colnames includes PS")
# }
# ExpressionSetInput <- ExpressionSet
# }
#
# # Make GeneID the row name if not already
# if("GeneID" %in% colnames(ExpressionSetInput)){
# ExpressionSetInput <-
# tibble::column_to_rownames(ExpressionSetInput, var = "GeneID")
# }
# # test
#
# m <- match.fun(boot_function)
#
# # Function to calculate mean
# average_function <- function(data, indices) {
# sampled_data <- data[indices, ]
# means <- apply(sampled_data, 2, m)
# return(means)
# }
#
# # Function to calculate standard deviation based on means
# std_dev_from_means <- function(boot_means) {
# return(apply(boot_means, 2, sd))
# }
#
# # Perform bootstrap resampling for means
# bootstrap_averages <- boot(data = ExpressionSetInput, statistic = average_function, R = permutations)
#
# # Calculate standard deviations from bootstrapped means
# bootstrap_std_dev <- std_dev_from_means(bootstrap_averages$t)
#
# return(bootstrap_std_dev)
# }I will just plot the relative expression. While it takes the mean
expression for a PS for each stage sum(Exp_{PS,Stage})/n
(where n is the number of genes in that PS) divided by the
maximum mean expression for a PS across stages,
i.e. (sum(Exp_{PS,Stage})/n)/(max(sum(Exp_{PS,Stage})/n)),
it is in fact equivalent to
(sum(Exp_{PS,Stage}))/(max(sum(Exp_{PS,Stage}))) (that is,
the relative sum expression level of genes in a PS).
myTAI::plotRE()Fs_PES %>%
myTAI::PlotRE(Groups = list(c(1:6),c(7,8)), main = "Fs; TPM")
Fd_PES %>%
myTAI::PlotRE(Groups = list(c(1:6),c(7,8)), main = "Fd; TPM")
Fs_PES.sqrt %>%
myTAI::PlotRE(Groups = list(c(1:6),c(7,8)), main = "Fs; sqrt(TPM)")
Fd_PES.sqrt %>%
myTAI::PlotRE(Groups = list(c(1:6),c(7,8)), main = "Fd; sqrt(TPM)")
Fs_PES.log2 %>%
myTAI::PlotRE(Groups = list(c(1:6),c(7,8)), main = "Fs; log2(TPM+1)")
Fd_PES.log2 %>%
myTAI::PlotRE(Groups = list(c(1:6),c(7,8)), main = "Fd; log2(TPM+1)")
Fs_PES.rank %>%
myTAI::PlotRE(Groups = list(c(1:6),c(7,8)), main = "Fs; rank(TPM)")
Fd_PES.rank %>%
myTAI::PlotRE(Groups = list(c(1:6),c(7,8)), main = "Fd; rank(TPM)")
Fs_PES.rlog %>%
myTAI::PlotRE(Groups = list(c(1:6),c(7,8)), main = "Fs; rlog(TPM)")
Fd_PES.rlog %>%
myTAI::PlotRE(Groups = list(c(1:6),c(7,8)), main = "Fd; rlog(TPM)")
# Fs_PES.log2 %>%
# tidyr::pivot_longer(!c(PS,GeneID), names_to = "Stage", values_to = "Expression") %>%
# dplyr::group_by(PS, Stage) %>%
# dplyr::summarise(Average = base::mean(Expression),
# y25 = quantile(Expression, 0.25),
# y75 = quantile(Expression, 0.75)) %>%
# ggplot2::ggplot(aes(x = Stage, y = Average, group = PS, colour = as.factor(PS))) +
# ggplot2::geom_point() +
# ggplot2::geom_errorbar(aes(ymin = y25, ymax = y75)) +
# ggplot2::facet_grid(. ~ PS) +
# ggplot2::theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1))Alternative function for more control.
plot_relative_expression <- function(PES, average = "mean", plot = T){
averaging_function <- base::match.fun(average)
PES_in <-
PES %>%
tidyr::pivot_longer(!c(PS, GeneID), values_to = "Expression", names_to = "Stage") %>%
dplyr::group_by(PS,Stage) %>%
dplyr::summarise(mean_expression = averaging_function(Expression),
n_genes = n())
max_expression <-
PES_in %>%
dplyr::group_by(PS) %>%
dplyr::summarise(max_expression = max(mean_expression),
min_expression = min(mean_expression))
# normalise expression between 0 and 1
PES_relative_expression <-
left_join(PES_in, max_expression, by = "PS") %>%
dplyr::mutate(relative_expression = (mean_expression-min_expression)/(max_expression-min_expression))
PES_relative_expression$Stage <-
factor(PES_relative_expression$Stage, levels = unique(colnames(PES)[-c(1,2)]))
# return if plotting is not desired.
if(!plot){return(PES_relative_expression)}
PES_relative_expression <-
dplyr::mutate(PES_relative_expression, PS_category = case_when(
PS < 7 ~ "ancient genes",
TRUE~ "young genes"
)) %>%
dplyr::mutate(PS_representativeness = case_when(
n_genes < 1000 ~ "low",
TRUE~ "high"
))
ggplot2::ggplot(PES_relative_expression,
aes(x = Stage,
y = relative_expression,
group = PS,
# linetype = PS_representativeness,
# linewidth = log2(n_genes),
colour = factor(PS))) +
# ggplot2::geom_line(alpha = 0.5) +
ggplot2::geom_line(linewidth = 2, linetype = "twodash") +
ggplot2::facet_wrap(vars(PS_category)) +
ggplot2::scale_color_brewer(palette="Dark2") +
ggplot2::labs(y = "relative expression", colour = "PS")
}plot_relative_expression(Fs_PES) + ggplot2::labs(title = "Fs", subtitle = "TPM")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_relative_expression(Fd_PES) + ggplot2::labs(title = "Fd", subtitle = "TPM")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_relative_expression(Fs_PES.sqrt) + ggplot2::labs(title = "Fs", subtitle = "sqrt(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_relative_expression(Fd_PES.sqrt) + ggplot2::labs(title = "Fd", subtitle = "sqrt(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_relative_expression(Fs_PES.log2) + ggplot2::labs(title = "Fs", subtitle = "log2(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_relative_expression(Fd_PES.log2) + ggplot2::labs(title = "Fd", subtitle = "log2(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_relative_expression(Fs_PES.rank) + ggplot2::labs(title = "Fs", subtitle = "rank(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_relative_expression(Fd_PES.rank) + ggplot2::labs(title = "Fd", subtitle = "rank(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_relative_expression(Fs_PES.rlog) + ggplot2::labs(title = "Fs", subtitle = "rlog(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_relative_expression(Fd_PES.rlog) + ggplot2::labs(title = "Fd", subtitle = "rlog(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
Thus, we see a general pattern that the hourglass pattern comes from the down-regulation of younger genes rather than an upregulation of older genes, especially for the sqrt and log2 transformed data. The dataset with ranked show that if we emphasise genes with low expression in the resulting TAI profile, we see increased expression of older genes in the mid-embryo. However, this transformation is non-standard and it is unclear what is even meant by “relative expression” in the context of the rank data.
Plot mean CI
Second attempt.
plot_PS_CImean <- function(PES){
PES_in <-
PES %>%
tidyr::pivot_longer(!c(PS, GeneID), values_to = "Expression", names_to = "Stage") %>%
dplyr::group_by(PS,Stage) %>%
dplyr::summarise(ggplot2::mean_cl_boot(Expression))
PES_in$Stage <-
factor(PES_in$Stage, levels = unique(colnames(PES)[-c(1,2)]))
PES_in %>%
ggplot2::ggplot(aes(x = Stage, y = y, colour = factor(PS), group = PS)) +
ggplot2::geom_ribbon(aes(ymin = ymin, ymax = ymax),
alpha = 0.1, color = NA) +
ggplot2::geom_line(size = 1) +
ggplot2::scale_color_brewer(palette="Dark2") +
ggplot2::facet_wrap(vars(PS), scales = "free_y") +
ggplot2::labs(y = "Mean expression level", colour = "PS")
}# Fd_PES %>%
# tidyr::pivot_longer(!c(PS, GeneID), values_to = "Expression", names_to = "Stage") %>%
# dplyr::group_by(PS,Stage) %>%
# dplyr::summarise(ggplot2::mean_cl_boot(Expression)) %>%
# ggplot2::ggplot(aes(x = Stage, y = y, colour = factor(PS), group = PS)) +
# ggplot2::geom_ribbon(aes(ymin = ymin, ymax = ymax),
# alpha = 0.1, color = NA) +
# ggplot2::geom_line(size = 1) +
# ggplot2::scale_color_brewer(palette="Dark2") +
# ggplot2::facet_wrap(vars(PS), scales = "free_y") +
# ggplot2::theme_classic()
plot_PS_CImean(Fs_PES) + ggplot2::labs(title = "Fs", subtitle = "TPM")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
plot_PS_CImean(Fd_PES) + ggplot2::labs(title = "Fd", subtitle = "TPM")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_PS_CImean(Fs_PES.sqrt) + ggplot2::labs(title = "Fs", subtitle = "sqrt(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_PS_CImean(Fd_PES.sqrt) + ggplot2::labs(title = "Fd", subtitle = "sqrt(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_PS_CImean(Fs_PES.log2) + ggplot2::labs(title = "Fs", subtitle = "log2(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_PS_CImean(Fd_PES.log2) + ggplot2::labs(title = "Fd", subtitle = "log2(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_PS_CImean(Fs_PES.rank) + ggplot2::labs(title = "Fs", subtitle = "rank(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_PS_CImean(Fd_PES.rank) + ggplot2::labs(title = "Fd", subtitle = "rank(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_PS_CImean(Fs_PES.rlog) + ggplot2::labs(title = "Fs", subtitle = "rlog(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
plot_PS_CImean(Fd_PES.rlog) + ggplot2::labs(title = "Fd", subtitle = "rlog(TPM)")## `summarise()` has grouped output by 'PS'. You can override using the `.groups`
## argument.
So what we see with more details (if the bootstrapping for the mean makes sense, which I do since the distribution underlying the expression is non-normal) is that the older genes are expressed with more variance and the youngest genes are completely downregulated in the ‘phylotypic stage’. Data is based on 95% confidence interval.
Get session info.
devtools::session_info()## ─ Session info ───────────────────────────────────────────────────────────────
## setting value
## version R version 4.2.2 (2022-10-31)
## os macOS Big Sur ... 10.16
## system x86_64, darwin17.0
## ui X11
## language (EN)
## collate en_US.UTF-8
## ctype en_US.UTF-8
## tz Europe/Berlin
## date 2024-02-14
## pandoc 3.1.6.2 @ /usr/local/bin/ (via rmarkdown)
##
## ─ Packages ───────────────────────────────────────────────────────────────────
## package * version date (UTC) lib source
## abind 1.4-5 2016-07-21 [1] CRAN (R 4.2.0)
## backports 1.4.1 2021-12-13 [1] CRAN (R 4.2.0)
## base64enc 0.1-3 2015-07-28 [1] CRAN (R 4.2.0)
## bit 4.0.5 2022-11-15 [1] CRAN (R 4.2.0)
## bit64 4.0.5 2020-08-30 [1] CRAN (R 4.2.0)
## broom 1.0.5 2023-06-09 [1] CRAN (R 4.2.0)
## bslib 0.5.1 2023-08-11 [1] CRAN (R 4.2.0)
## cachem 1.0.8 2023-05-01 [1] CRAN (R 4.2.0)
## callr 3.7.3 2022-11-02 [1] CRAN (R 4.2.0)
## car 3.1-2 2023-03-30 [1] CRAN (R 4.2.0)
## carData 3.0-5 2022-01-06 [1] CRAN (R 4.2.0)
## checkmate 2.2.0 2023-04-27 [1] CRAN (R 4.2.0)
## cli 3.6.1 2023-03-23 [1] CRAN (R 4.2.0)
## cluster 2.1.4 2022-08-22 [1] CRAN (R 4.2.2)
## codetools 0.2-19 2023-02-01 [1] CRAN (R 4.2.0)
## colorspace 2.1-0 2023-01-23 [1] CRAN (R 4.2.0)
## cowplot 1.1.1 2020-12-30 [1] CRAN (R 4.2.0)
## crayon 1.5.2 2022-09-29 [1] CRAN (R 4.2.0)
## data.table 1.14.8 2023-02-17 [1] CRAN (R 4.2.0)
## devtools 2.4.5 2022-10-11 [1] CRAN (R 4.2.0)
## digest 0.6.33 2023-07-07 [1] CRAN (R 4.2.0)
## dplyr * 1.1.3 2023-09-03 [1] CRAN (R 4.2.0)
## ellipsis 0.3.2 2021-04-29 [1] CRAN (R 4.2.0)
## evaluate 0.22 2023-09-29 [1] CRAN (R 4.2.2)
## fansi 1.0.5 2023-10-08 [1] CRAN (R 4.2.2)
## farver 2.1.1 2022-07-06 [1] CRAN (R 4.2.0)
## fastmap 1.1.1 2023-02-24 [1] CRAN (R 4.2.0)
## forcats * 1.0.0 2023-01-29 [1] CRAN (R 4.2.0)
## foreach 1.5.2 2022-02-02 [1] CRAN (R 4.2.0)
## foreign 0.8-85 2023-09-09 [1] CRAN (R 4.2.0)
## Formula 1.2-5 2023-02-24 [1] CRAN (R 4.2.0)
## fs 1.6.3 2023-07-20 [1] CRAN (R 4.2.0)
## generics 0.1.3 2022-07-05 [1] CRAN (R 4.2.0)
## ggplot2 * 3.4.4 2023-10-12 [1] CRAN (R 4.2.2)
## ggpubr 0.6.0 2023-02-10 [1] CRAN (R 4.2.0)
## ggsci 3.0.0 2023-03-08 [1] CRAN (R 4.2.0)
## ggsignif 0.6.4 2022-10-13 [1] CRAN (R 4.2.0)
## glue 1.6.2 2022-02-24 [1] CRAN (R 4.2.0)
## gridExtra 2.3 2017-09-09 [1] CRAN (R 4.2.0)
## gtable 0.3.4 2023-08-21 [1] CRAN (R 4.2.0)
## Hmisc 5.1-1 2023-09-12 [1] CRAN (R 4.2.0)
## hms 1.1.3 2023-03-21 [1] CRAN (R 4.2.0)
## htmlTable 2.4.1 2022-07-07 [1] CRAN (R 4.2.0)
## htmltools 0.5.6.1 2023-10-06 [1] CRAN (R 4.2.2)
## htmlwidgets 1.6.2 2023-03-17 [1] CRAN (R 4.2.0)
## httpuv 1.6.11 2023-05-11 [1] CRAN (R 4.2.2)
## iterators 1.0.14 2022-02-05 [1] CRAN (R 4.2.0)
## jquerylib 0.1.4 2021-04-26 [1] CRAN (R 4.2.0)
## jsonlite 1.8.7 2023-06-29 [1] CRAN (R 4.2.0)
## knitr 1.44 2023-09-11 [1] CRAN (R 4.2.0)
## labeling 0.4.3 2023-08-29 [1] CRAN (R 4.2.0)
## later 1.3.1 2023-05-02 [1] CRAN (R 4.2.2)
## lattice 0.21-9 2023-10-01 [1] CRAN (R 4.2.2)
## lifecycle 1.0.3 2022-10-07 [1] CRAN (R 4.2.0)
## lubridate * 1.9.3 2023-09-27 [1] CRAN (R 4.2.0)
## magrittr 2.0.3 2022-03-30 [1] CRAN (R 4.2.0)
## Matrix 1.5-4.1 2023-05-18 [1] CRAN (R 4.2.0)
## memoise 2.0.1 2021-11-26 [1] CRAN (R 4.2.0)
## mime 0.12 2021-09-28 [1] CRAN (R 4.2.0)
## miniUI 0.1.1.1 2018-05-18 [1] CRAN (R 4.2.0)
## munsell 0.5.0 2018-06-12 [1] CRAN (R 4.2.0)
## myTAI * 1.0.1.9000 2023-12-07 [1] Github (drostlab/myTAI@27a2639)
## nnet 7.3-19 2023-05-03 [1] CRAN (R 4.2.2)
## pillar 1.9.0 2023-03-22 [1] CRAN (R 4.2.0)
## pkgbuild 1.4.2 2023-06-26 [1] CRAN (R 4.2.0)
## pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 4.2.0)
## pkgload 1.3.3 2023-09-22 [1] CRAN (R 4.2.0)
## plyr 1.8.9 2023-10-02 [1] CRAN (R 4.2.2)
## prettyunits 1.2.0 2023-09-24 [1] CRAN (R 4.2.0)
## processx 3.8.2 2023-06-30 [1] CRAN (R 4.2.0)
## profvis 0.3.8 2023-05-02 [1] CRAN (R 4.2.0)
## promises 1.2.1 2023-08-10 [1] CRAN (R 4.2.2)
## ps 1.7.5 2023-04-18 [1] CRAN (R 4.2.0)
## purrr * 1.0.2 2023-08-10 [1] CRAN (R 4.2.2)
## R6 2.5.1 2021-08-19 [1] CRAN (R 4.2.0)
## RColorBrewer 1.1-3 2022-04-03 [1] CRAN (R 4.2.0)
## Rcpp 1.0.11 2023-07-06 [1] CRAN (R 4.2.0)
## readr * 2.1.4 2023-02-10 [1] CRAN (R 4.2.0)
## remotes 2.4.2.1 2023-07-18 [1] CRAN (R 4.2.2)
## reshape2 1.4.4 2020-04-09 [1] CRAN (R 4.2.0)
## rlang 1.1.1 2023-04-28 [1] CRAN (R 4.2.0)
## rmarkdown 2.25 2023-09-18 [1] CRAN (R 4.2.2)
## rpart 4.1.21 2023-10-09 [1] CRAN (R 4.2.2)
## rstatix 0.7.2 2023-02-01 [1] CRAN (R 4.2.0)
## rstudioapi 0.15.0 2023-07-07 [1] CRAN (R 4.2.0)
## sass 0.4.7 2023-07-15 [1] CRAN (R 4.2.0)
## scales 1.2.1 2022-08-20 [1] CRAN (R 4.2.0)
## sessioninfo 1.2.2 2021-12-06 [1] CRAN (R 4.2.0)
## shiny 1.7.5.1 2023-10-14 [1] CRAN (R 4.2.2)
## stringi 1.7.12 2023-01-11 [1] CRAN (R 4.2.0)
## stringr * 1.5.0 2022-12-02 [1] CRAN (R 4.2.0)
## tibble * 3.2.1 2023-03-20 [1] CRAN (R 4.2.0)
## tidyr * 1.3.0 2023-01-24 [1] CRAN (R 4.2.0)
## tidyselect 1.2.0 2022-10-10 [1] CRAN (R 4.2.0)
## tidyverse * 2.0.0 2023-02-22 [1] CRAN (R 4.2.0)
## timechange 0.2.0 2023-01-11 [1] CRAN (R 4.2.0)
## tzdb 0.4.0 2023-05-12 [1] CRAN (R 4.2.2)
## urlchecker 1.0.1 2021-11-30 [1] CRAN (R 4.2.0)
## usethis 2.2.2 2023-07-06 [1] CRAN (R 4.2.0)
## utf8 1.2.3 2023-01-31 [1] CRAN (R 4.2.0)
## vctrs 0.6.4 2023-10-12 [1] CRAN (R 4.2.2)
## vroom 1.6.4 2023-10-02 [1] CRAN (R 4.2.2)
## withr 2.5.1 2023-09-26 [1] CRAN (R 4.2.0)
## xfun 0.40 2023-08-09 [1] CRAN (R 4.2.2)
## xtable 1.8-4 2019-04-21 [1] CRAN (R 4.2.0)
## yaml 2.3.7 2023-01-23 [1] CRAN (R 4.2.0)
##
## [1] /Library/Frameworks/R.framework/Versions/4.2/Resources/library
##
## ──────────────────────────────────────────────────────────────────────────────